/*


!
!  Dalton, a molecular electronic structure program
!  Copyright (C) The Dalton Authors (see AUTHORS file for details).
!
!  This program is free software; you can redistribute it and/or
!  modify it under the terms of the GNU Lesser General Public
!  License version 2.1 as published by the Free Software Foundation.
!
!  This program is distributed in the hope that it will be useful,
!  but WITHOUT ANY WARRANTY; without even the implied warranty of
!  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
!  Lesser General Public License for more details.
!
!  If a copy of the GNU LGPL v2.1 was not distributed with this
!  code, you can obtain one at https://www.gnu.org/licenses/old-licenses/lgpl-2.1.en.html.
!

!

*/
/*-*-mode: C; c-indentation-style: "bsd"; c-basic-offset: 4; -*-*/
/* fun-p86c.c:

   REFERENCE: J.P. Perdew, Phys. Rev. B, 33, 8822 (1986).

   Perdew local (P86) correlation functional, ie VWN with gradient correction.

   Automatically generated code implementing P86c functional and its
   derivatives. Generated by func-codegen.pl being a part of a 
   "Automatic code generation framework for analytical functional
    derivative evaluation", Pawel Salek, 2004
    
    This functional is connected by making following changes:
    1. add "extern Functional p86cFunctional;" to 'functionals.h'
    2. add "&p86cFunctional," to 'functionals.c'
    3. add "fun-p86c}.c" to 'Makefile.am', 'Makefile.in' or 'Makefile'.

    This functional has been generated from following input:
    ------ cut here -------
c0:    1.667e-3;
cn1:   2.568e-3;
cn2:   2.3266e-2;
cn3:   7.389e-6;
cd1:   1.0;
cd2:   8.723;
cd3:   0.472;
cd4:   1.0e4*cn3;
cinf:  c0+cn1/cd1;
phi1:  (9.0*%PI)^(1/6); 
phi2:  0.11;

rho:  rhoa + rhob;
grad: sqrt(grada*grada + gradb*gradb + 2*gradab);
zeta: (rhoa-rhob)/(rhoa+rhob);
dzet: 2^(1/3)*sqrt(((1+zeta)/2)^(5/3) + ((1-zeta)/2)^(5/3));

rs:    (3/(4*%PI*rho))^(1/3);
crho:  c0 + (cn1+cn2*rs+cn3*rs^2)/(1+cd2*rs+cd3*(rs^2)+cd4*(rs^3));
phi:   phi1*phi2*(cinf/crho)*grad*(rho^(-7.0/6.0));

K(rhoa,rhob,grada,gradb,gradab):=exp(-phi)*crho*(grad^2)*(rho^(-4.0/3.0))/dzet;

    ------ cut here -------
*/

 
/* strictly conform to XOPEN ANSI C standard */
#if !defined(SYS_DEC)
/* XOPEN compliance is missing on old Tru64 4.0E Alphas and pow() prototype
 * is not specified. */
#define _XOPEN_SOURCE          500
#define _XOPEN_SOURCE_EXTENDED 1
#endif
#include <math.h>
#include <stddef.h>
#include "general.h"

#define __CVERSION__
 
#include "functionals.h"
 
/* INTERFACE PART */
static integer p86c_isgga(void) { return 1; } /* FIXME: detect! */
static integer p86c_read(const char *conf_line);
static real p86c_energy(const FunDensProp *dp);
static void p86c_first(FunFirstFuncDrv *ds,   real factor,
                         const FunDensProp *dp);
static void p86c_second(FunSecondFuncDrv *ds, real factor,
                          const FunDensProp *dp);
static void p86c_third(FunThirdFuncDrv *ds,   real factor,
                         const FunDensProp *dp);
 
Functional P86cFunctional = {
  "P86c",       /* name */
  p86c_isgga,   /* gga-corrected */
   1,
  p86c_read,
  NULL,
  p86c_energy,
  p86c_first,
  p86c_second,
  p86c_third
};
 
/* IMPLEMENTATION PART */
static integer
p86c_read(const char *conf_line)
{
    fun_set_hf_weight(0);
    return 1;
}

static real
p86c_energy(const FunDensProp *dp)
{
    real res;
    real rhoa = dp->rhoa, rhob = dp->rhob;
    real grada = dp->grada, gradb = dp->gradb, gradab = dp->gradab;

    real t1, t2, t3, t4, t5, t6, t7, t8, t9, t10;
    real t11, t12, t13, t14;

    t1 = pow(gradb,2.0)+2.0*gradab+pow(grada,2.0);
    t2 = rhob+rhoa;
    t3 = 1/pow(2.0,0.66666666666667);
    t4 = rhoa-1.0*rhob;
    t5 = 1/t2;
    t6 = pow(3.0,0.66666666666667);
    t7 = 1/pow(4.0,0.66666666666667);
    t8 = 1/pow(3.141592653589793,0.66666666666667);
    t9 = 1/pow(t2,0.66666666666667);
    t10 = pow(3.0,0.33333333333333);
    t11 = 1/pow(4.0,0.33333333333333);
    t12 = 1/pow(3.141592653589793,0.33333333333333);
    t13 = 1/pow(t2,0.33333333333333);
    t14 = (0.023266*t10*t11*t12*t13+7.389E-6*t6*t7*t8*t9+
        0.002568)/(0.472*t6*t7*t8*t9+0.01763993811759*t5+8.723000000000001*
        t10*t11*t12*t13+1.0)+0.001667;

   /* code */
    res = t1*t14/(pow(2.0,0.33333333333333)*pow(2.718281828459045,
        6.718719623277062E-4*pow(3.141592653589793,0.16666666666667)*
        sqrt(t1)/(t14*pow(t2,1.166666666666667)))*pow(t2,1.333333333333333)*
        sqrt(0.5*t3*pow(t4*t5+1.0,1.666666666666667)+0.5*t3*pow(1.0-
        1.0*t4*t5,1.666666666666667)));

    return res;
}

static void
p86c_first(FunFirstFuncDrv *ds, real factor, const FunDensProp *dp)
{
    real dfdra, dfdrb, dfdga, dfdgb, dfdgab;
    real rhoa = dp->rhoa, rhob = dp->rhob;
    real grada = dp->grada, gradb = dp->gradb, gradab = dp->gradab;

    real t1, t2, t3, t4, t5, t6, t7, t8, t9, t10;
    real t11, t12, t13, t14, t15, t16, t17, t18;
    real t19, t20, t21, t22, t23, t24, t25, t26;
    real t27, t28, t29, t30, t31, t32, t33, t34;
    real t35, t36, t37, t38, t39, t40, t41, t42;
    real t43;

    t1 = 1/pow(2.0,0.33333333333333);
    t2 = pow(gradb,2.0)+2.0*gradab+pow(grada,2.0);
    t3 = rhob+rhoa;
    t4 = 1/pow(t3,1.333333333333333);
    t5 = 1/pow(2.0,0.66666666666667);
    t6 = rhoa-1.0*rhob;
    t7 = 1/t3;
    t8 = 1.0-1.0*t6*t7;
    t9 = t6*t7+1.0;
    t10 = sqrt(0.5*t5*pow(t9,1.666666666666667)+0.5*t5*pow(t8,
        1.666666666666667));
    t11 = 1/t10;
    t12 = 0.31830988618379;
    t13 = 1/pow(t3,2.0);
    t14 = pow(3.0,0.66666666666667);
    t15 = 1/pow(4.0,0.66666666666667);
    t16 = 1/pow(3.141592653589793,0.66666666666667);
    t17 = 1/pow(t3,1.666666666666667);
    t18 = pow(3.0,0.33333333333333);
    t19 = 1/pow(4.0,0.33333333333333);
    t20 = 1/pow(3.141592653589793,0.33333333333333);
    t21 = 1/pow(t3,1.333333333333333);
    t22 = 1/pow(t3,0.66666666666667);
    t23 = 1/pow(t3,0.33333333333333);
    t24 = 0.023266*t18*t19*t20*t23+7.389E-6*t14*t15*t16*t22+
        0.002568;
    t25 = 0.0554175*t12*t7+8.723000000000001*t18*t19*t20*
        t23+0.472*t14*t15*t16*t22+1.0;
    t26 = 1/t25;
    t27 = (-0.00775533333333*t18*t19*t20*t21-4.926E-6*t14*
        t15*t16*t17)*t26-1.0*(-2.907666666666667*t18*t19*t20*t21-0.31466666666667*
        t14*t15*t16*t17-0.0554175*t12*t13)*t24/pow(t25,2.0);
    t28 = pow(3.141592653589793,0.16666666666667);
    t29 = sqrt(t2);
    t30 = 1/pow(t3,1.166666666666667);
    t31 = t24*t26+0.001667;
    t32 = 1/t31;
    t33 = 1/pow(2.718281828459045,6.718719623277062E-4*t28*
        t29*t30*t32);
    t34 = t1*t2*t4*t11*t27*t33;
    t35 = t6*t13;
    t36 = -1.0*t7;
    t37 = pow(t8,0.66666666666667);
    t38 = -1.0*t13*t6;
    t39 = pow(t9,0.66666666666667);
    t40 = 1/pow(t10,3.0);
    t41 = -1.333333333333333*t1*t11*t2*t31*t33/pow(t3,2.333333333333333);
    t42 = t1*
        t11*t2*t31*(7.838506227156574E-4*t28*t29*t32/pow(t3,2.166666666666667)+
        6.718719623277062E-4*t27*t28*t29*t30/pow(t31,2.0))*t33*t4;
    t43 = 1/
        pow(t3,2.5);

   /* code */
    dfdra = -0.5*t1*t2*t31*t33*t4*t40*(0.83333333333333*t39*
        t5*(t7+t38)+0.83333333333333*(t36+t35)*t37*t5)+t42+t41+t34;
    dfdrb = -
        0.5*t1*t2*t31*t33*t4*t40*(0.83333333333333*t37*t5*(t7+t35)+
        0.83333333333333*(t36+t38)*t39*t5)+t42+t41+t34;
    dfdga = 2.0*t1*t11*t31*t33*t4*grada-6.718719623277062E-4*
        t1*t28*grada*t29*t43*t11*t33;
    dfdgb = 2.0*t1*t11*t31*t33*t4*gradb-6.718719623277062E-4*
        t1*t28*gradb*t29*t43*t11*t33;
    dfdgab = 2.0*t1*t11*t31*t33*t4-6.718719623277062E-4*t1*
        t28*t29*t43*t11*t33;


   /* Final assignment */
    ds->df1000 += factor*dfdra;
    ds->df0100 += factor*dfdrb;
    ds->df0010 += factor*dfdga;
    ds->df0001 += factor*dfdgb;
    ds->df00001+= factor*dfdgab;
}

static void
p86c_second(FunSecondFuncDrv *ds, real factor, const FunDensProp *dp)
{
    real dfdra, dfdrb, dfdga, dfdgb, dfdgab;
    real d2fdrara, d2fdrarb, d2fdraga, d2fdragb, d2fdrbrb, d2fdraab, 
         d2fdrbab, d2fdgaga, d2fdgbgb, d2fdrbga, d2fdrbgb,
         d2fdragab,d2fdrbgab, d2fdgagb, d2fdgagab, d2fdgbgab, d2fdgabgab;
    real rhoa = dp->rhoa, rhob = dp->rhob;
    real grada = dp->grada, gradb = dp->gradb, gradab = dp->gradab;

    real t1, t2, t3, t4, t5, t6, t7, t8, t9, t10;
    real t11, t12, t13, t14, t15, t16, t17, t18;
    real t19, t20, t21, t22, t23, t24, t25, t26;
    real t27, t28, t29, t30, t31, t32, t33, t34;
    real t35, t36, t37, t38, t39, t40, t41, t42;
    real t43, t44, t45, t46, t47, t48, t49, t50;
    real t51, t52, t53, t54, t55, t56, t57, t58;
    real t59, t60, t61, t62, t63, t64, t65, t66;
    real t67, t68, t69, t70, t71, t72, t73, t74;
    real t75, t76, t77, t78, t79, t80, t81, t82;
    real t83, t84, t85, t86, t87, t88, t89, t90;
    real t91, t92, t93, t94, t95, t96, t97, t98;
    real t99, t100, t101, t102, t103, t104, t105;

    t1 = 1/pow(2.0,0.33333333333333);
    t2 = pow(grada,2.0);
    t3 = pow(gradb,2.0);
    t4 = 2.0*gradab+t3+t2;
    t5 = rhob+rhoa;
    t6 = 1/pow(t5,1.333333333333333);
    t7 = 1/pow(2.0,0.66666666666667);
    t8 = rhoa-1.0*rhob;
    t9 = 1/t5;
    t10 = 1.0-1.0*t8*t9;
    t11 = t8*t9+1.0;
    t12 = sqrt(0.5*pow(t11,1.666666666666667)*t7+0.5*pow(t10,
        1.666666666666667)*t7);
    t13 = 1/t12;
    t14 = 0.31830988618379;
    t15 = 1/pow(t5,2.0);
    t16 = pow(3.0,0.66666666666667);
    t17 = 1/pow(4.0,0.66666666666667);
    t18 = 1/pow(3.141592653589793,0.66666666666667);
    t19 = 1/pow(t5,1.666666666666667);
    t20 = pow(3.0,0.33333333333333);
    t21 = 1/pow(4.0,0.33333333333333);
    t22 = pow(3.141592653589793,0.33333333333333);
    t23 = 1/t22;
    t24 = 1/pow(t5,1.333333333333333);
    t25 = -2.907666666666667*t20*t21*t23*t24-0.31466666666667*
        t16*t17*t18*t19-0.0554175*t14*t15;
    t26 = 1/pow(t5,0.66666666666667);
    t27 = 1/pow(t5,0.33333333333333);
    t28 = 0.023266*t20*t21*t23*t27+7.389E-6*t16*t17*t18*t26+
        0.002568;
    t29 = 0.0554175*t14*t9+8.723000000000001*t20*t21*t23*
        t27+0.472*t16*t17*t18*t26+1.0;
    t30 = 1/pow(t29,2.0);
    t31 = -0.00775533333333*t20*t21*t23*t24-4.926E-6*t16*
        t17*t18*t19;
    t32 = 1/t29;
    t33 = t31*t32-1.0*t25*t28*t30;
    t34 = pow(3.141592653589793,0.16666666666667);
    t35 = sqrt(t4);
    t36 = 1/pow(t5,1.166666666666667);
    t37 = t28*t32+0.001667;
    t38 = 1/t37;
    t39 = 1/pow(2.718281828459045,6.718719623277062E-4*t34*
        t35*t36*t38);
    t40 = t1*t4*t6*t13*t33*t39;
    t41 = t8*t15;
    t42 = -1.0*t9;
    t43 = t42+t41;
    t44 = pow(t10,0.66666666666667);
    t45 = -1.0*t15*t8;
    t46 = t9+t45;
    t47 = pow(t11,0.66666666666667);
    t48 = 0.83333333333333*t46*t47*t7+0.83333333333333*t43*
        t44*t7;
    t49 = 1/pow(t12,3.0);
    t50 = 1/pow(t5,2.333333333333333);
    t51 = -1.333333333333333*t1*t4*t50*t13*t37*t39;
    t52 = 1/pow(t37,2.0);
    t53 = 1/pow(t5,2.166666666666667);
    t54 = 7.838506227156574E-4*t34*t35*t53*t38+6.718719623277062E-4*
        t34*t35*t36*t33*t52;
    t55 = t1*t4*t6*t13*t37*t54*t39;
    t56 = t9+t41;
    t57 = t42+t45;
    t58 = 0.83333333333333*t47*t57*t7+0.83333333333333*t44*
        t56*t7;
    t59 = 1/pow(t5,2.5);
    t60 = -6.718719623277062E-4*t1*t34*t35*t59*t13*t39;
    t61 = 2.0*
        t1*t13*t37*t39*t6;
    t62 = 1/pow(t5,3.0);
    t63 = 1/pow(t5,2.666666666666667);
    t64 = 1/pow(t5,2.333333333333334);
    t65 = -1.0*t28*t30*(3.876888888888889*t20*t21*t23*t64+
        0.52444444444444*t16*t17*t18*t63+0.110835*t14*t62)+(0.01034044444444*
        t20*t21*t23*t64+8.210000000000001E-6*t16*t17*t18*t63)*t32-
        2.0*t25*t30*t31+2.0*pow(t25,2.0)*t28/pow(t29,3.0);
    t66 = t1*t4*t6*t13*t65*t39;
    t67 = -2.666666666666667*t1*t4*t50*t13*t33*t39;
    t68 = 1/pow(t12,5.0);
    t69 = 1/pow(t10,0.33333333333333);
    t70 = -2.0*t62*t8;
    t71 = 2.0*t15;
    t72 = 1/pow(t11,0.33333333333333);
    t73 = 2.0*t62*t8;
    t74 = -2.0*t15;
    t75 = 3.111111111111111*t1*t13*t37*t39*t4/pow(t5,3.333333333333333);
    t76 = t1*
        t13*t37*t39*t4*(6.718719623277062E-4*t34*t35*t36*t65*t52-0.00156770124543*
        t34*t35*t53*t33*t52-0.00169834301588*t34*t35*t38/pow(t5,3.166666666666667)-
        0.00134374392466*pow(t33,2.0)*t34*t35*t36/pow(t37,3.0))*t6;
    t77 = 2.0*
        t1*t13*t33*t39*t4*t54*t6;
    t78 = -2.666666666666667*t1*t4*t50*t13*t37*t54*t39;
    t79 = t1*
        t13*t37*t39*t4*pow(t54,2.0)*t6;
    t80 = 1/pow(t5,3.5);
    t81 = 8.958292831036083E-4*t1*t34*grada*t35*t80*t13*t39;
    t82 = 2.0*
        t1*t13*t33*t39*t6*grada;
    t83 = -6.718719623277062E-4*t1*t34*grada*t35*t59*t13*
        t33*t38*t39;
    t84 = -2.666666666666667*t1*grada*t50*t13*t37*t39;
    t85 = 1/t35;
    t86 = t1*t4*t6*t13*t37*(7.838506227156574E-4*t34*grada*
        t85*t53*t38+6.718719623277062E-4*t34*grada*t85*t36*t33*t52)*
        t39;
    t87 = -6.718719623277062E-4*t1*t34*grada*t35*t59*t13*
        t54*t39;
    t88 = 2.0*t1*t13*t37*t39*t54*t6*grada;
    t89 = 8.958292831036083E-4*t1*t34*gradb*t35*t80*t13*t39;
    t90 = 2.0*
        t1*t13*t33*t39*t6*gradb;
    t91 = -6.718719623277062E-4*t1*t34*gradb*t35*t59*t13*
        t33*t38*t39;
    t92 = -2.666666666666667*t1*gradb*t50*t13*t37*t39;
    t93 = t1*t4*t6*t13*t37*(7.838506227156574E-4*t34*gradb*
        t85*t53*t38+6.718719623277062E-4*t34*gradb*t85*t36*t33*t52)*
        t39;
    t94 = -6.718719623277062E-4*t1*t34*gradb*t35*t59*t13*
        t54*t39;
    t95 = 2.0*t1*t13*t37*t39*t54*t6*gradb;
    t96 = 8.958292831036083E-4*t1*t34*t35*t80*t13*t39;
    t97 = 2.0*t1*t13*t33*t39*t6;
    t98 = -6.718719623277062E-4*t1*t34*t35*t59*t13*t33*t38*
        t39;
    t99 = -2.666666666666667*t1*t50*t13*t37*t39;
    t100 = t1*t4*t6*t13*t37*(7.838506227156574E-4*t34*t85*
        t53*t38+6.718719623277062E-4*t34*t85*t36*t33*t52)*t39;
    t101 = -6.718719623277062E-4*t1*t34*t35*t59*t13*t54*t39;
    t102 = 2.0*
        t1*t13*t37*t39*t54*t6;
    t103 = t99+t98+t97+t96-1.0*t1*t37*t39*t48*t49*t6+3.359359811638531E-4*
        t1*t34*t35*t59*t48*t49*t39+t102+t101+t100;
    t104 = t99+t98+t97+t96-1.0*t1*t37*t39*t49*t58*t6+3.359359811638531E-4*
        t1*t34*t35*t59*t58*t49*t39+t102+t101+t100;
    t105 = 1/pow(t5,3.666666666666667);

   /* code */
    dfdra = -0.5*t1*t37*t39*t4*t48*t49*t6+t55+t51+t40;
    dfdrb = -0.5*t1*t37*t39*t4*t49*t58*t6+t55+t51+t40;
    dfdga = 2.0*t1*t13*t37*t39*t6*grada-6.718719623277062E-4*
        t1*t34*grada*t35*t59*t13*t39;
    dfdgb = 2.0*t1*t13*t37*t39*t6*gradb-6.718719623277062E-4*
        t1*t34*gradb*t35*t59*t13*t39;
    dfdgab = t61+t60;
    d2fdrara = t79+t78+t77+t76+t75-0.5*t1*t37*t39*t4*t49*
        t6*(0.83333333333333*t47*t7*(t74+t73)+0.55555555555556*pow(t46,
        2.0)*t7*t72+0.83333333333333*t44*t7*(t71+t70)+0.55555555555556*
        pow(t43,2.0)*t69*t7)+0.75*t1*t37*t39*t4*pow(t48,2.0)*t6*t68+
        t67+t66-1.0*t1*t37*t39*t4*t48*t49*t54*t6-1.0*t1*t33*t39*t4*
        t48*t49*t6+1.333333333333333*t1*t4*t50*t48*t49*t37*t39;
    d2fdrarb = -0.5*t1*t37*t39*t4*t49*t6*(1.666666666666667*
        t47*t62*t7*t8-1.666666666666667*t44*t62*t7*t8+0.55555555555556*
        t46*t57*t7*t72+0.55555555555556*t43*t56*t69*t7)+t79+t78+t77+
        t76+t75+0.75*t1*t37*t39*t4*t48*t58*t6*t68+t67+t66-0.5*t1*t37*
        t39*t4*t49*t54*t58*t6-0.5*t1*t33*t39*t4*t49*t58*t6-0.5*t1*
        t37*t39*t4*t48*t49*t54*t6-0.5*t1*t33*t39*t4*t48*t49*t6+0.66666666666667*
        t1*t4*t50*t58*t49*t37*t39+0.66666666666667*t1*t4*t50*t48*t49*
        t37*t39;
    d2fdraga = -1.0*t1*t37*t39*t48*t49*t6*grada+t88+t87+t86+
        t84+t83+t82+t81+3.359359811638531E-4*t1*t34*grada*t35*t59*
        t48*t49*t39;
    d2fdragb = -1.0*t1*t37*t39*t48*t49*t6*gradb+t95+t94+t93+
        t92+t91+t90+t89+3.359359811638531E-4*t1*t34*gradb*t35*t59*
        t48*t49*t39;
    d2fdrbrb = t79+t78+t77+t76+t75-0.5*t1*t37*t39*t4*t49*
        t6*(0.83333333333333*t44*t7*(t74+t70)+0.55555555555556*pow(t57,
        2.0)*t7*t72+0.83333333333333*t47*t7*(t71+t73)+0.55555555555556*
        pow(t56,2.0)*t69*t7)+0.75*t1*t37*t39*t4*pow(t58,2.0)*t6*t68+
        t67+t66-1.0*t1*t37*t39*t4*t49*t54*t58*t6-1.0*t1*t33*t39*t4*
        t49*t58*t6+1.333333333333333*t1*t4*t50*t58*t49*t37*t39;
    d2fdraab = t103;
    d2fdrbab = t104;
    d2fdgaga = t61+4.514119337620827E-7*t1*t22*t2*t105*t13*
        t38*t39+t60-0.00201561588698*t1*t34*t2*t85*t59*t13*t39;
    d2fdgbgb = t61+4.514119337620827E-7*t1*t22*t3*t105*t13*
        t38*t39+t60-0.00201561588698*t1*t34*t3*t85*t59*t13*t39;
    d2fdrbga = -1.0*t1*t37*t39*t49*t58*t6*grada+t88+t87+t86+
        t84+t83+t82+t81+3.359359811638531E-4*t1*t34*grada*t35*t59*
        t58*t49*t39;
    d2fdrbgb = -1.0*t1*t37*t39*t49*t58*t6*gradb+t95+t94+t93+
        t92+t91+t90+t89+3.359359811638531E-4*t1*t34*gradb*t35*t59*
        t58*t49*t39;
    d2fdgagb = 4.514119337620827E-7*t1*t22*grada*gradb*t105*
        t13*t38*t39-0.00201561588698*t1*t34*grada*gradb*t85*t59*t13*
        t39;
    d2fdragab = t103;
    d2fdrbgab = t104;
    d2fdgagab = 4.514119337620827E-7*t1*t22*grada*t105*t13*
        t38*t39-0.00201561588698*t1*t34*grada*t85*t59*t13*t39;
    d2fdgbgab = 4.514119337620827E-7*t1*t22*gradb*t105*t13*
        t38*t39-0.00201561588698*t1*t34*gradb*t85*t59*t13*t39;
    d2fdgabgab = 4.514119337620827E-7*t1*t22*t105*t13*t38*
        t39-0.00201561588698*t1*t34*t85*t59*t13*t39;


    ds->df1000 += factor*dfdra;
    ds->df0100 += factor*dfdrb;
    ds->df0010 += factor*dfdga;
    ds->df0001 += factor*dfdgb;
    ds->df00001+= factor*dfdgab;

    ds->df2000 += factor*d2fdrara;
    ds->df1100 += factor*d2fdrarb;
    ds->df1010 += factor*d2fdraga;
    ds->df1001 += factor*d2fdragb;
    ds->df10001+= factor*d2fdraab;
    ds->df0200 += factor*d2fdrbrb;
    ds->df0110 += factor*d2fdrbga;
    ds->df0101 += factor*d2fdrbgb;
    ds->df01001+= factor*d2fdrbgab;
    ds->df0020 += factor*d2fdgaga;
    ds->df0011 += factor*d2fdgagb;
    ds->df00101+= factor*d2fdgagab;
    ds->df0002 += factor*d2fdgbgb;
    ds->df00011+= factor*d2fdgbgab;
    ds->df00002+= factor*d2fdgabgab;
}

static void
p86c_third(FunThirdFuncDrv *ds, real factor, const FunDensProp *dp)
{
    real dfdra, dfdrb, dfdga, dfdgb, dfdgab;
    real d2fdrara, d2fdrarb, d2fdraga, d2fdragb, d2fdrbrb, d2fdraab, 
         d2fdrbab, d2fdgaga, d2fdgbgb, d2fdrbga, d2fdrbgb,
         d2fdragab,d2fdrbgab, d2fdgagb, d2fdgagab, d2fdgbgab, d2fdgabgab;
    real d3fdrararb, d3fdraraga, d3fdraragb, d3fdrbrbab, d3fdraraab,
         d3fdrarbrb, d3fdrarbga, d3fdrarbgb, d3fdrarbab, d3fdragaga,
         d3fdragbgb, d3fdrarara, d3fdrbrbrb, d3fdrbrbga, d3fdrbrbgb,
         d3fdrbgaga, d3fdrbgbga, d3fdrbgbgb, d3fdgagaga,
         d3fdragagb, d3fdrbgagb, d3fdgagagb, d3fdgagbgb, d3fdgbgbgb;
    real rhoa = dp->rhoa, rhob = dp->rhob;
    real grada = dp->grada, gradb = dp->gradb, gradab = dp->gradab;

    real t1, t2, t3, t4, t5, t6, t7, t8, t9, t10;
    real t11, t12, t13, t14, t15, t16, t17, t18;
    real t19, t20, t21, t22, t23, t24, t25, t26;
    real t27, t28, t29, t30, t31, t32, t33, t34;
    real t35, t36, t37, t38, t39, t40, t41, t42;
    real t43, t44, t45, t46, t47, t48, t49, t50;
    real t51, t52, t53, t54, t55, t56, t57, t58;
    real t59, t60, t61, t62, t63, t64, t65, t66;
    real t67, t68, t69, t70, t71, t72, t73, t74;
    real t75, t76, t77, t78, t79, t80, t81, t82;
    real t83, t84, t85, t86, t87, t88, t89, t90;
    real t91, t92, t93, t94, t95, t96, t97, t98;
    real t99, t100, t101, t102, t103, t104, t105;
    real t106, t107, t108, t109, t110, t111, t112;
    real t113, t114, t115, t116, t117, t118, t119;
    real t120, t121, t122, t123, t124, t125, t126;
    real t127, t128, t129, t130, t131, t132, t133;
    real t134, t135, t136, t137, t138, t139, t140;
    real t141, t142, t143, t144, t145, t146, t147;
    real t148, t149, t150, t151, t152, t153, t154;
    real t155, t156, t157, t158, t159, t160, t161;
    real t162, t163, t164, t165, t166, t167, t168;
    real t169, t170, t171, t172, t173, t174, t175;
    real t176, t177, t178, t179, t180, t181, t182;
    real t183, t184, t185, t186, t187, t188, t189;
    real t190, t191, t192, t193, t194, t195, t196;
    real t197, t198, t199, t200, t201, t202, t203;
    real t204, t205, t206, t207, t208, t209, t210;
    real t211, t212, t213, t214, t215, t216, t217;
    real t218, t219, t220, t221, t222, t223, t224;
    real t225, t226, t227, t228, t229, t230, t231;
    real t232, t233, t234, t235, t236, t237, t238;
    real t239, t240, t241, t242, t243, t244, t245;
    real t246, t247, t248, t249, t250, t251, t252;
    real t253, t254, t255, t256, t257, t258, t259;
    real t260, t261, t262, t263, t264;

    t1 = 1/pow(2.0,0.33333333333333);
    t2 = pow(grada,2.0);
    t3 = pow(gradb,2.0);
    t4 = 2.0*gradab+t3+t2;
    t5 = rhob+rhoa;
    t6 = 1/pow(t5,1.333333333333333);
    t7 = 1/pow(2.0,0.66666666666667);
    t8 = rhoa-1.0*rhob;
    t9 = 1/t5;
    t10 = 1.0-1.0*t8*t9;
    t11 = t8*t9+1.0;
    t12 = sqrt(0.5*pow(t11,1.666666666666667)*t7+0.5*pow(t10,
        1.666666666666667)*t7);
    t13 = 1/t12;
    t14 = 0.31830988618379;
    t15 = 1/pow(t5,2.0);
    t16 = pow(3.0,0.66666666666667);
    t17 = 1/pow(4.0,0.66666666666667);
    t18 = 1/pow(3.141592653589793,0.66666666666667);
    t19 = 1/pow(t5,1.666666666666667);
    t20 = pow(3.0,0.33333333333333);
    t21 = 1/pow(4.0,0.33333333333333);
    t22 = pow(3.141592653589793,0.33333333333333);
    t23 = 1/t22;
    t24 = 1/pow(t5,1.333333333333333);
    t25 = -2.907666666666667*t20*t21*t23*t24-0.31466666666667*
        t16*t17*t18*t19-0.0554175*t14*t15;
    t26 = 1/pow(t5,0.66666666666667);
    t27 = 1/pow(t5,0.33333333333333);
    t28 = 0.023266*t20*t21*t23*t27+7.389E-6*t16*t17*t18*t26+
        0.002568;
    t29 = 0.0554175*t14*t9+8.723000000000001*t20*t21*t23*
        t27+0.472*t16*t17*t18*t26+1.0;
    t30 = 1/pow(t29,2.0);
    t31 = -0.00775533333333*t20*t21*t23*t24-4.926E-6*t16*
        t17*t18*t19;
    t32 = 1/t29;
    t33 = t31*t32-1.0*t25*t28*t30;
    t34 = pow(3.141592653589793,0.16666666666667);
    t35 = sqrt(t4);
    t36 = 1/pow(t5,1.166666666666667);
    t37 = t28*t32+0.001667;
    t38 = 1/t37;
    t39 = 1/pow(2.718281828459045,6.718719623277062E-4*t34*
        t35*t36*t38);
    t40 = t1*t4*t6*t13*t33*t39;
    t41 = t8*t15;
    t42 = -1.0*t9;
    t43 = t42+t41;
    t44 = pow(t10,0.66666666666667);
    t45 = -1.0*t15*t8;
    t46 = t9+t45;
    t47 = pow(t11,0.66666666666667);
    t48 = 0.83333333333333*t46*t47*t7+0.83333333333333*t43*
        t44*t7;
    t49 = 1/pow(t12,3.0);
    t50 = 1/pow(t5,2.333333333333333);
    t51 = -1.333333333333333*t1*t4*t50*t13*t37*t39;
    t52 = 1/pow(t37,2.0);
    t53 = 1/pow(t5,2.166666666666667);
    t54 = 7.838506227156574E-4*t34*t35*t53*t38+6.718719623277062E-4*
        t34*t35*t36*t33*t52;
    t55 = t1*t4*t6*t13*t37*t54*t39;
    t56 = t9+t41;
    t57 = t42+t45;
    t58 = 0.83333333333333*t47*t57*t7+0.83333333333333*t44*
        t56*t7;
    t59 = 1/pow(t5,2.5);
    t60 = -6.718719623277062E-4*t1*t34*t35*t59*t13*t39;
    t61 = 2.0*
        t1*t13*t37*t39*t6;
    t62 = pow(t25,2.0);
    t63 = 1/pow(t29,3.0);
    t64 = 1/pow(t5,3.0);
    t65 = 1/pow(t5,2.666666666666667);
    t66 = 1/pow(t5,2.333333333333334);
    t67 = 3.876888888888889*t20*t21*t23*t66+0.52444444444444*
        t16*t17*t18*t65+0.110835*t14*t64;
    t68 = 0.01034044444444*t20*t21*t23*t66+8.210000000000001E-6*
        t16*t17*t18*t65;
    t69 = -1.0*t28*t30*t67+2.0*t28*t62*t63+t68*t32-2.0*t25*
        t30*t31;
    t70 = t1*t4*t6*t13*t69*t39;
    t71 = -2.666666666666667*t1*t4*t50*t13*t33*t39;
    t72 = pow(t48,2.0);
    t73 = 1/pow(t12,5.0);
    t74 = pow(t43,2.0);
    t75 = 1/pow(t10,0.33333333333333);
    t76 = -2.0*t64*t8;
    t77 = 2.0*t15;
    t78 = t77+t76;
    t79 = pow(t46,2.0);
    t80 = 1/pow(t11,0.33333333333333);
    t81 = 2.0*t64*t8;
    t82 = -2.0*t15;
    t83 = t82+t81;
    t84 = 0.83333333333333*t47*t7*t83+0.55555555555556*t7*
        t79*t80+0.83333333333333*t44*t7*t78+0.55555555555556*t7*t74*
        t75;
    t85 = 1/pow(t5,3.333333333333333);
    t86 = 3.111111111111111*t1*t4*t85*t13*t37*t39;
    t87 = pow(t33,2.0);
    t88 = 1/pow(t37,3.0);
    t89 = 1/pow(t5,3.166666666666667);
    t90 = -0.00169834301588*t34*t35*t89*t38-0.00156770124543*
        t34*t35*t53*t33*t52+6.718719623277062E-4*t34*t35*t36*t69*t52-
        0.00134374392466*t34*t35*t36*t87*t88;
    t91 = t1*t4*t6*t13*t37*t90*t39;
    t92 = 2.0*t1*t13*t33*t39*t4*t54*t6;
    t93 = -2.666666666666667*t1*t4*t50*t13*t37*t54*t39;
    t94 = pow(t54,
        2.0);
    t95 = t1*t4*t6*t13*t37*t94*t39;
    t96 = 0.55555555555556*t46*t57*t7*t80+1.666666666666667*
        t47*t64*t7*t8-1.666666666666667*t44*t64*t7*t8+0.55555555555556*
        t43*t56*t7*t75;
    t97 = 1/pow(t5,3.5);
    t98 = 8.958292831036083E-4*t1*t34*grada*t35*t97*t13*t39;
    t99 = 2.0*
        t1*t13*t33*t39*t6*grada;
    t100 = -6.718719623277062E-4*t1*t34*grada*t35*t59*t13*
        t33*t38*t39;
    t101 = -2.666666666666667*t1*grada*t50*t13*t37*t39;
    t102 = 1/
        t35;
    t103 = 7.838506227156574E-4*t34*grada*t102*t53*t38+6.718719623277062E-4*
        t34*grada*t102*t36*t33*t52;
    t104 = t1*t4*t6*t13*t37*t103*t39;
    t105 = -6.718719623277062E-4*t1*t34*grada*t35*t59*t13*
        t54*t39;
    t106 = 2.0*t1*t13*t37*t39*t54*t6*grada;
    t107 = 8.958292831036083E-4*t1*t34*gradb*t35*t97*t13*
        t39;
    t108 = 2.0*t1*t13*t33*t39*t6*gradb;
    t109 = -6.718719623277062E-4*t1*t34*gradb*t35*t59*t13*
        t33*t38*t39;
    t110 = -2.666666666666667*t1*gradb*t50*t13*t37*t39;
    t111 = 7.838506227156574E-4*t34*gradb*t102*t53*t38+6.718719623277062E-4*
        t34*gradb*t102*t36*t33*t52;
    t112 = t1*t4*t6*t13*t37*t111*t39;
    t113 = -6.718719623277062E-4*t1*t34*gradb*t35*t59*t13*
        t54*t39;
    t114 = 2.0*t1*t13*t37*t39*t54*t6*gradb;
    t115 = pow(t58,2.0);
    t116 = pow(t56,2.0);
    t117 = t82+t76;
    t118 = pow(t57,2.0);
    t119 = t77+t81;
    t120 = 0.55555555555556*t118*t7*t80+0.55555555555556*
        t116*t7*t75+0.83333333333333*t119*t47*t7+0.83333333333333*
        t117*t44*t7;
    t121 = 3.359359811638531E-4*t1*t34*t35*t59*t48*t49*t39;
    t122 = 8.958292831036083E-4*t1*t34*t35*t97*t13*t39;
    t123 = 2.0*
        t1*t13*t33*t39*t6;
    t124 = -6.718719623277062E-4*t1*t34*t35*t59*t13*t33*t38*
        t39;
    t125 = -1.0*t1*t37*t39*t48*t49*t6;
    t126 = -2.666666666666667*t1*t50*t13*t37*t39;
    t127 = 6.718719623277062E-4*t34*t102*t36*t33*t52;
    t128 = 7.838506227156574E-4*t34*t102*t53*t38;
    t129 = t128+t127;
    t130 = t1*t4*t6*t13*t37*t129*t39;
    t131 = -6.718719623277062E-4*t1*t34*t35*t59*t13*t54*t39;
    t132 = 2.0*
        t1*t13*t37*t39*t54*t6;
    t133 = t132+t131+t130+t126+t125+t124+t123+t122+t121;
    t134 = 3.359359811638531E-4*t1*t34*t35*t59*t58*t49*t39;
    t135 = -
        1.0*t1*t37*t39*t49*t58*t6;
    t136 = t132+t131+t130+t126+t135+t124+t123+t122+t134;
    t137 = 1/
        pow(t5,3.666666666666667);
    t138 = -0.00201561588698*t1*t34*grada*t102*t59*t13*t39;
    t139 = 4.514119337620827E-7*t1*t22*grada*t137*t13*t38*
        t39;
    t140 = -0.00201561588698*t1*t34*gradb*t102*t59*t13*t39;
    t141 = 4.514119337620827E-7*t1*t22*gradb*t137*t13*t38*
        t39;
    t142 = 1/pow(t5,4.0);
    t143 = 1/pow(t5,3.666666666666667);
    t144 = 1/pow(t5,3.333333333333334);
    t145 = -3.0*t25*t30*t68+6.0*t25*t28*t63*t67-3.0*t30*t31*
        t67+6.0*t31*t62*t63+(-0.0241277037037*t20*t21*t23*t144-2.189333333333334E-5*
        t16*t17*t18*t143)*t32-1.0*(-9.046074074074074*t20*t21*t23*
        t144-1.398518518518519*t16*t17*t18*t143-0.332505*t14*t142)*
        t28*t30-6.0*pow(t25,3.0)*t28/pow(t29,4.0);
    t146 = t1*t4*t6*t13*t145*t39;
    t147 = -4.0*t1*t4*t50*t13*t69*t39;
    t148 = 1.5*t1*t33*t39*t4*t48*t58*t6*t73;
    t149 = -1.0*t1*t33*t39*t4*t49*t6*t96;
    t150 = 9.333333333333332*t1*t4*t85*t13*t33*t39;
    t151 = 1/pow(t12,7.0);
    t152 = -2.0*t1*t4*t50*t58*t48*t73*t37*t39;
    t153 = 1.333333333333333*t1*t4*t50*t96*t49*t37*t39;
    t154 = 1/
        pow(t10,1.333333333333333);
    t155 = 6.0*t142*t8;
    t156 = 1/pow(t11,1.333333333333333);
    t157 = -6.0*t142*t8;
    t158 = -10.37037037037037*t1*t13*t37*t39*t4/pow(t5,4.333333333333333);
    t159 = t1*
        t13*t37*t39*t4*t6*(0.00470310373629*t34*t35*t53*t87*t88-0.00403123177397*
        t34*t35*t36*t69*t33*t88-0.00235155186815*t34*t35*t53*t69*t52+
        0.00509502904765*t34*t35*t89*t33*t52+6.718719623277062E-4*
        t34*t35*t36*t145*t52+0.00537808621697*t34*t35*t38/pow(t5,4.166666666666667)+
        0.00403123177397*pow(t33,3.0)*t34*t35*t36/pow(t37,4.0));
    t160 = 3.0*
        t1*t13*t33*t39*t4*t6*t90;
    t161 = -4.0*t1*t4*t50*t13*t37*t90*t39;
    t162 = 3.0*t1*t13*t39*t4*t54*t6*t69;
    t163 = -8.0*t1*t4*t50*t13*t33*t54*t39;
    t164 = 1.5*t1*t37*t39*t4*t48*t54*t58*t6*t73;
    t165 = -1.0*t1*t37*t39*t4*t49*t54*t6*t96;
    t166 = 9.333333333333332*t1*t4*t85*t13*t37*t54*t39;
    t167 = 3.0*
        t1*t13*t37*t39*t4*t54*t6*t90;
    t168 = 3.0*t1*t13*t33*t39*t4*t6*t94;
    t169 = -4.0*t1*t4*t50*t13*t37*t94*t39;
    t170 = t1*t13*t37*t39*t4*pow(t54,3.0)*t6;
    t171 = 1/pow(t5,4.5);
    t172 = -0.00209026832724*t1*t34*grada*t35*t171*t13*t39;
    t173 = 2.0*
        t1*t13*t39*t6*t69*grada;
    t174 = -5.333333333333333*t1*grada*t50*t13*t33*t39;
    t175 = -
        6.718719623277062E-4*t1*t34*grada*t35*t59*t13*t69*t38*t39;
    t176 = 0.00179165856621*
        t1*t34*grada*t35*t97*t13*t33*t38*t39;
    t177 = 6.222222222222221*t1*grada*t85*t13*t37*t39;
    t178 = t1*t4*t6*t13*t37*(-0.00169834301588*t34*grada*
        t102*t89*t38-0.00156770124543*t34*grada*t102*t53*t33*t52+6.718719623277062E-4*
        t34*grada*t102*t36*t69*t52-0.00134374392466*t34*grada*t102*
        t36*t87*t88)*t39;
    t179 = -6.718719623277062E-4*t1*t34*grada*t35*t59*t13*
        t90*t39;
    t180 = 2.0*t1*t13*t37*t39*t6*t90*grada;
    t181 = 2.0*t1*t103*t13*t33*t39*t4*t6;
    t182 = -2.666666666666667*t1*t4*t50*t13*t37*t103*t39;
    t183 = 0.00179165856621*
        t1*t34*grada*t35*t97*t13*t54*t39;
    t184 = 4.0*t1*t13*t33*t39*t54*t6*grada;
    t185 = -0.00134374392466*t1*t34*grada*t35*t59*t13*t33*
        t38*t54*t39;
    t186 = -5.333333333333333*t1*grada*t50*t13*t37*t54*t39;
    t187 = 2.0*
        t1*t103*t13*t37*t39*t4*t54*t6;
    t188 = -6.718719623277062E-4*t1*t34*grada*t35*t59*t13*
        t94*t39;
    t189 = 2.0*t1*t13*t37*t39*t6*t94*grada;
    t190 = -0.00209026832724*t1*t34*gradb*t35*t171*t13*t39;
    t191 = 2.0*
        t1*t13*t39*t6*t69*gradb;
    t192 = -5.333333333333333*t1*gradb*t50*t13*t33*t39;
    t193 = -
        6.718719623277062E-4*t1*t34*gradb*t35*t59*t13*t69*t38*t39;
    t194 = 0.00179165856621*
        t1*t34*gradb*t35*t97*t13*t33*t38*t39;
    t195 = 6.222222222222221*t1*gradb*t85*t13*t37*t39;
    t196 = t1*t4*t6*t13*t37*(-0.00169834301588*t34*gradb*
        t102*t89*t38-0.00156770124543*t34*gradb*t102*t53*t33*t52+6.718719623277062E-4*
        t34*gradb*t102*t36*t69*t52-0.00134374392466*t34*gradb*t102*
        t36*t87*t88)*t39;
    t197 = -6.718719623277062E-4*t1*t34*gradb*t35*t59*t13*
        t90*t39;
    t198 = 2.0*t1*t13*t37*t39*t6*t90*gradb;
    t199 = 2.0*t1*t111*t13*t33*t39*t4*t6;
    t200 = -2.666666666666667*t1*t4*t50*t13*t37*t111*t39;
    t201 = 0.00179165856621*
        t1*t34*gradb*t35*t97*t13*t54*t39;
    t202 = 4.0*t1*t13*t33*t39*t54*t6*gradb;
    t203 = -0.00134374392466*t1*t34*gradb*t35*t59*t13*t33*
        t38*t54*t39;
    t204 = -5.333333333333333*t1*gradb*t50*t13*t37*t54*t39;
    t205 = 2.0*
        t1*t111*t13*t37*t39*t4*t54*t6;
    t206 = -6.718719623277062E-4*t1*t34*gradb*t35*t59*t13*
        t94*t39;
    t207 = 2.0*t1*t13*t37*t39*t6*t94*gradb;
    t208 = -0.00209026832724*t1*t34*t35*t171*t13*t39;
    t209 = 2.0*t1*t13*t39*t6*t69;
    t210 = -5.333333333333333*t1*t50*t13*t33*t39;
    t211 = -6.718719623277062E-4*t1*t34*t35*t59*t13*t69*t38*
        t39;
    t212 = 0.00179165856621*t1*t34*t35*t97*t13*t33*t38*t39;
    t213 = 6.222222222222221*
        t1*t85*t13*t37*t39;
    t214 = t1*t4*t6*t13*t37*(-0.00169834301588*t34*t102*t89*
        t38-0.00156770124543*t34*t102*t53*t33*t52+6.718719623277062E-4*
        t34*t102*t36*t69*t52-0.00134374392466*t34*t102*t36*t87*t88)*
        t39;
    t215 = -6.718719623277062E-4*t1*t34*t35*t59*t13*t90*t39;
    t216 = 2.0*
        t1*t13*t37*t39*t6*t90;
    t217 = 2.0*t1*t129*t13*t33*t39*t4*t6;
    t218 = -2.666666666666667*t1*t4*t50*t13*t37*t129*t39;
    t219 = 0.00179165856621*
        t1*t34*t35*t97*t13*t54*t39;
    t220 = 4.0*t1*t13*t33*t39*t54*t6;
    t221 = -0.00134374392466*t1*t34*t35*t59*t13*t33*t38*t54*
        t39;
    t222 = -5.333333333333333*t1*t50*t13*t37*t54*t39;
    t223 = 2.0*t1*t129*t13*t37*t39*t4*t54*t6;
    t224 = -6.718719623277062E-4*t1*t34*t35*t59*t13*t94*t39;
    t225 = 2.0*
        t1*t13*t37*t39*t6*t94;
    t226 = 0.00268748784931*t1*t34*t2*t102*t97*t13*t39;
    t227 = 4.514119337620827E-7*t1*t22*t2*t137*t13*t33*t52*
        t39;
    t228 = 1/pow(t5,4.666666666666667);
    t229 = -6.018825783494436E-7*t1*t22*t2*t228*t13*t38*t39;
    t230 = -
        0.00201561588698*t1*t34*t2*t102*t59*t13*t33*t38*t39;
    t231 = 1/pow(t35,3.0);
    t232 = t1*t4*t6*t13*t37*(t128-7.838506227156574E-4*t34*
        t2*t231*t53*t38+t127-6.718719623277062E-4*t34*t2*t231*t36*
        t33*t52)*t39;
    t233 = -0.00134374392466*t1*t34*grada*t35*t59*t13*t103*
        t39;
    t234 = 4.0*t1*t103*t13*t37*t39*t6*grada;
    t235 = -0.00201561588698*t1*t34*t2*t102*t59*t13*t54*t39;
    t236 = 4.514119337620827E-7*t1*t22*t2*t137*t13*t38*t54*
        t39;
    t237 = 0.00268748784931*t1*t34*grada*gradb*t102*t97*t13*
        t39;
    t238 = 4.514119337620827E-7*t1*t22*grada*gradb*t137*t13*
        t33*t52*t39;
    t239 = -6.018825783494436E-7*t1*t22*grada*gradb*t228*
        t13*t38*t39;
    t240 = -0.00201561588698*t1*t34*grada*gradb*t102*t59*
        t13*t33*t38*t39;
    t241 = t1*t4*t6*t13*t37*(-7.838506227156574E-4*t34*grada*
        gradb*t231*t53*t38-6.718719623277062E-4*t34*grada*gradb*t231*
        t36*t33*t52)*t39;
    t242 = -6.718719623277062E-4*t1*t34*gradb*t35*t59*t13*
        t103*t39;
    t243 = 2.0*t1*t103*t13*t37*t39*t6*gradb;
    t244 = -6.718719623277062E-4*t1*t34*grada*t35*t59*t13*
        t111*t39;
    t245 = 2.0*t1*t111*t13*t37*t39*t6*grada;
    t246 = -0.00201561588698*t1*t34*grada*gradb*t102*t59*
        t13*t54*t39;
    t247 = 4.514119337620827E-7*t1*t22*grada*gradb*t137*t13*
        t38*t54*t39;
    t248 = 0.00268748784931*t1*t34*t3*t102*t97*t13*t39;
    t249 = 4.514119337620827E-7*t1*t22*t3*t137*t13*t33*t52*
        t39;
    t250 = -6.018825783494436E-7*t1*t22*t3*t228*t13*t38*t39;
    t251 = -
        0.00201561588698*t1*t34*t3*t102*t59*t13*t33*t38*t39;
    t252 = t1*t4*t6*t13*t37*(t128-7.838506227156574E-4*t34*
        t3*t231*t53*t38+t127-6.718719623277062E-4*t34*t3*t231*t36*
        t33*t52)*t39;
    t253 = -0.00134374392466*t1*t34*gradb*t35*t59*t13*t111*
        t39;
    t254 = 4.0*t1*t111*t13*t37*t39*t6*gradb;
    t255 = -0.00201561588698*t1*t34*t3*t102*t59*t13*t54*t39;
    t256 = 4.514119337620827E-7*t1*t22*t3*t137*t13*t38*t54*
        t39;
    t257 = -6.0*t64;
    t258 = 6.0*t64;
    t259 = t132+t236+t131+t235+t234+t233+t232+t126+t135+t124+
        t230+t229-2.257059668810414E-7*t1*t22*t2*t137*t58*t49*t38*
        t39+t227+t123+t122+t226+t134+0.00100780794349*t1*t34*t2*t102*
        t59*t58*t49*t39;
    t260 = pow(grada,3.0);
    t261 = 1.772453850905516;
    t262 = 1/pow(t5,4.833333333333334);
    t263 = 1/t4;
    t264 = pow(gradb,3.0);

   /* code */
    dfdra = -0.5*t1*t37*t39*t4*t48*t49*t6+t55+t51+t40;
    dfdrb = -0.5*t1*t37*t39*t4*t49*t58*t6+t55+t51+t40;
    dfdga = 2.0*t1*t13*t37*t39*t6*grada-6.718719623277062E-4*
        t1*t34*grada*t35*t59*t13*t39;
    dfdgb = 2.0*t1*t13*t37*t39*t6*gradb-6.718719623277062E-4*
        t1*t34*gradb*t35*t59*t13*t39;
    dfdgab = t61+t60;
    d2fdrara = t95+t93+t92+t91+t86-0.5*t1*t37*t39*t4*t49*
        t6*t84+0.75*t1*t37*t39*t4*t6*t72*t73+t71+t70-1.0*t1*t37*t39*
        t4*t48*t49*t54*t6-1.0*t1*t33*t39*t4*t48*t49*t6+1.333333333333333*
        t1*t4*t50*t48*t49*t37*t39;
    d2fdrarb = -0.5*t1*t37*t39*t4*t49*t6*t96+t95+t93+t92+
        t91+t86+0.75*t1*t37*t39*t4*t48*t58*t6*t73+t71+t70-0.5*t1*t37*
        t39*t4*t49*t54*t58*t6-0.5*t1*t33*t39*t4*t49*t58*t6-0.5*t1*
        t37*t39*t4*t48*t49*t54*t6-0.5*t1*t33*t39*t4*t48*t49*t6+0.66666666666667*
        t1*t4*t50*t58*t49*t37*t39+0.66666666666667*t1*t4*t50*t48*t49*
        t37*t39;
    d2fdraga = -1.0*t1*t37*t39*t48*t49*t6*grada+t99+t98+3.359359811638531E-4*
        t1*t34*grada*t35*t59*t48*t49*t39+t106+t105+t104+t101+t100;
    d2fdragb = -
        1.0*t1*t37*t39*t48*t49*t6*gradb+3.359359811638531E-4*t1*t34*
        gradb*t35*t59*t48*t49*t39+t114+t113+t112+t110+t109+t108+t107;
    d2fdrbrb = t95+
        t93+t92+t91+t86+0.75*t1*t115*t37*t39*t4*t6*t73+t71+t70-1.0*
        t1*t37*t39*t4*t49*t54*t58*t6-1.0*t1*t33*t39*t4*t49*t58*t6-
        0.5*t1*t120*t37*t39*t4*t49*t6+1.333333333333333*t1*t4*t50*
        t58*t49*t37*t39;
    d2fdraab = t133;
    d2fdrbab = t136;
    d2fdgaga = t61+4.514119337620827E-7*t1*t22*t2*t137*t13*
        t38*t39+t60-0.00201561588698*t1*t34*t2*t102*t59*t13*t39;
    d2fdgbgb = t61+
        4.514119337620827E-7*t1*t22*t3*t137*t13*t38*t39+t60-0.00201561588698*
        t1*t34*t3*t102*t59*t13*t39;
    d2fdrbga = -1.0*t1*t37*t39*t49*t58*t6*grada+t99+t98+3.359359811638531E-4*
        t1*t34*grada*t35*t59*t58*t49*t39+t106+t105+t104+t101+t100;
    d2fdrbgb = -
        1.0*t1*t37*t39*t49*t58*t6*gradb+3.359359811638531E-4*t1*t34*
        gradb*t35*t59*t58*t49*t39+t114+t113+t112+t110+t109+t108+t107;
    d2fdgagb = 4.514119337620827E-7*t1*t22*grada*gradb*t137*
        t13*t38*t39-0.00201561588698*t1*t34*grada*gradb*t102*t59*t13*
        t39;
    d2fdragab = t133;
    d2fdrbgab = t136;
    d2fdgagab = t139+t138;
    d2fdgbgab = t141+t140;
    d2fdgabgab = 4.514119337620827E-7*t1*t22*t137*t13*t38*
        t39-0.00201561588698*t1*t34*t102*t59*t13*t39;
    d3fdrararb = 1.5*t1*t37*t39*t4*t48*t6*t73*t96-0.5*t1*
        t37*t39*t4*t49*t58*t6*t94-1.0*t1*t37*t39*t4*t48*t49*t6*t94-
        0.5*t1*t37*t39*t4*t49*t58*t6*t90-1.0*t1*t37*t39*t4*t48*t49*
        t6*t90+0.75*t1*t37*t39*t4*t58*t6*t73*t84-0.5*t1*t37*t39*t4*
        t49*t54*t6*t84-0.5*t1*t33*t39*t4*t49*t6*t84-0.5*t1*t37*t39*
        t4*t49*t6*(0.55555555555556*t57*t7*t80*t83+2.222222222222222*
        t46*t64*t7*t8*t80-2.222222222222222*t43*t64*t7*t75*t8-0.18518518518519*
        t156*t57*t7*t79+0.55555555555556*t56*t7*t75*t78-0.18518518518519*
        t154*t56*t7*t74+0.83333333333333*t47*(2.0*t64+t157)*t7+0.83333333333333*
        t44*(t155-2.0*t64)*t7)+0.75*t1*t37*t39*t4*t54*t6*t72*t73+0.75*
        t1*t33*t39*t4*t6*t72*t73-1.875*t1*t151*t37*t39*t4*t58*t6*t72-
        0.5*t1*t39*t4*t49*t58*t6*t69-1.0*t1*t39*t4*t48*t49*t6*t69-
        1.0*t1*t33*t39*t4*t49*t54*t58*t6-2.0*t1*t33*t39*t4*t48*t49*
        t54*t6+1.333333333333333*t1*t4*t50*t58*t49*t37*t54*t39+2.666666666666667*
        t1*t4*t50*t48*t49*t37*t54*t39-1.0*t1*t4*t50*t72*t73*t37*t39+
        0.66666666666667*t1*t4*t50*t84*t49*t37*t39-1.555555555555555*
        t1*t4*t85*t58*t49*t37*t39-3.111111111111111*t1*t4*t85*t48*
        t49*t37*t39+1.333333333333333*t1*t4*t50*t58*t49*t33*t39+2.666666666666667*
        t1*t4*t50*t48*t49*t33*t39+t170+t169+t168+t167+t166+t165+t164+
        t163+t162+t161+t160+t159+t158+t153+t152+t150+t149+t148+t147+
        t146;
    d3fdraraga = -1.0*t1*t37*t39*t49*t6*t84*grada+1.5*t1*
        t37*t39*t6*t72*t73*grada-2.0*t1*t37*t39*t48*t49*t54*t6*grada-
        2.0*t1*t33*t39*t48*t49*t6*grada-1.0*t1*t103*t37*t39*t4*t48*
        t49*t6-5.039039717457796E-4*t1*t34*grada*t35*t59*t72*t73*t39+
        6.718719623277062E-4*t1*t34*grada*t35*t59*t48*t49*t54*t39+
        3.359359811638531E-4*t1*t34*grada*t35*t59*t84*t49*t39-8.958292831036083E-4*
        t1*t34*grada*t35*t97*t48*t49*t39+6.718719623277062E-4*t1*t34*
        grada*t35*t59*t48*t49*t33*t38*t39+2.666666666666667*t1*grada*
        t50*t48*t49*t37*t39+t189+t188+t187+t186+t185+t184+t183+t182+
        t181+t180+t179+t178+t177+t176+t175+t174+t173+t172;
    d3fdraragb = -1.0*t1*t37*t39*t49*t6*t84*gradb+1.5*t1*
        t37*t39*t6*t72*t73*gradb-2.0*t1*t37*t39*t48*t49*t54*t6*gradb-
        2.0*t1*t33*t39*t48*t49*t6*gradb-1.0*t1*t111*t37*t39*t4*t48*
        t49*t6-5.039039717457796E-4*t1*t34*gradb*t35*t59*t72*t73*t39+
        6.718719623277062E-4*t1*t34*gradb*t35*t59*t48*t49*t54*t39+
        3.359359811638531E-4*t1*t34*gradb*t35*t59*t84*t49*t39-8.958292831036083E-4*
        t1*t34*gradb*t35*t97*t48*t49*t39+6.718719623277062E-4*t1*t34*
        gradb*t35*t59*t48*t49*t33*t38*t39+2.666666666666667*t1*gradb*
        t50*t48*t49*t37*t39+t207+t206+t205+t204+t203+t202+t201+t200+
        t199+t198+t197+t196+t195+t194+t193+t192+t191+t190;
    d3fdrbrbab = 1.5*t1*t115*t37*t39*t6*t73-2.0*t1*t37*t39*
        t49*t54*t58*t6-1.0*t1*t129*t37*t39*t4*t49*t58*t6-2.0*t1*t33*
        t39*t49*t58*t6-1.0*t1*t120*t37*t39*t49*t6-5.039039717457796E-4*
        t1*t34*t35*t59*t115*t73*t39+6.718719623277062E-4*t1*t34*t35*
        t59*t58*t49*t54*t39-8.958292831036083E-4*t1*t34*t35*t97*t58*
        t49*t39+3.359359811638531E-4*t1*t34*t35*t59*t120*t49*t39+6.718719623277062E-4*
        t1*t34*t35*t59*t58*t49*t33*t38*t39+2.666666666666667*t1*t50*
        t58*t49*t37*t39+t225+t224+t223+t222+t221+t220+t219+t218+t217+
        t216+t215+t214+t213+t212+t211+t210+t209+t208;
    d3fdraraab = -1.0*t1*t37*t39*t49*t6*t84+1.5*t1*t37*t39*
        t6*t72*t73-2.0*t1*t37*t39*t48*t49*t54*t6-1.0*t1*t129*t37*t39*
        t4*t48*t49*t6-2.0*t1*t33*t39*t48*t49*t6-5.039039717457796E-4*
        t1*t34*t35*t59*t72*t73*t39+6.718719623277062E-4*t1*t34*t35*
        t59*t48*t49*t54*t39+3.359359811638531E-4*t1*t34*t35*t59*t84*
        t49*t39-8.958292831036083E-4*t1*t34*t35*t97*t48*t49*t39+6.718719623277062E-4*
        t1*t34*t35*t59*t48*t49*t33*t38*t39+2.666666666666667*t1*t50*
        t48*t49*t37*t39+t225+t224+t223+t222+t221+t220+t219+t218+t217+
        t216+t215+t214+t213+t212+t211+t210+t209+t208;
    d3fdrarbrb = 1.5*t1*t37*t39*t4*t58*t6*t73*t96-1.0*t1*
        t37*t39*t4*t49*t58*t6*t94-0.5*t1*t37*t39*t4*t48*t49*t6*t94-
        1.0*t1*t37*t39*t4*t49*t58*t6*t90-0.5*t1*t37*t39*t4*t48*t49*
        t6*t90-0.5*t1*t37*t39*t4*t49*t6*(2.222222222222222*t57*t64*
        t7*t8*t80+0.55555555555556*t119*t46*t7*t80-2.222222222222222*
        t56*t64*t7*t75*t8-5.0*t142*t47*t7*t8+5.0*t142*t44*t7*t8+0.55555555555556*
        t117*t43*t7*t75-1.666666666666667*t47*t64*t7+1.666666666666667*
        t44*t64*t7-0.18518518518519*t118*t156*t46*t7-0.18518518518519*
        t116*t154*t43*t7)+0.75*t1*t115*t37*t39*t4*t54*t6*t73+0.75*
        t1*t120*t37*t39*t4*t48*t6*t73+0.75*t1*t115*t33*t39*t4*t6*t73-
        1.0*t1*t39*t4*t49*t58*t6*t69-0.5*t1*t39*t4*t48*t49*t6*t69-
        2.0*t1*t33*t39*t4*t49*t54*t58*t6-1.0*t1*t33*t39*t4*t48*t49*
        t54*t6-0.5*t1*t120*t37*t39*t4*t49*t54*t6-0.5*t1*t120*t33*t39*
        t4*t49*t6-1.875*t1*t115*t151*t37*t39*t4*t48*t6+2.666666666666667*
        t1*t4*t50*t58*t49*t37*t54*t39+1.333333333333333*t1*t4*t50*
        t48*t49*t37*t54*t39-1.0*t1*t4*t50*t115*t73*t37*t39-3.111111111111111*
        t1*t4*t85*t58*t49*t37*t39-1.555555555555555*t1*t4*t85*t48*
        t49*t37*t39+0.66666666666667*t1*t4*t50*t120*t49*t37*t39+2.666666666666667*
        t1*t4*t50*t58*t49*t33*t39+1.333333333333333*t1*t4*t50*t48*
        t49*t33*t39+t170+t169+t168+t167+t166+t165+t164+t163+t162+t161+
        t160+t159+t158+t153+t152+t150+t149+t148+t147+t146;
    d3fdrarbga = -1.0*t1*t37*t39*t49*t6*t96*grada+1.5*t1*
        t37*t39*t48*t58*t6*t73*grada-1.0*t1*t37*t39*t49*t54*t58*t6*
        grada-1.0*t1*t33*t39*t49*t58*t6*grada-1.0*t1*t37*t39*t48*t49*
        t54*t6*grada-1.0*t1*t33*t39*t48*t49*t6*grada-0.5*t1*t103*t37*
        t39*t4*t49*t58*t6-0.5*t1*t103*t37*t39*t4*t48*t49*t6-5.039039717457796E-4*
        t1*t34*grada*t35*t59*t58*t48*t73*t39+3.359359811638531E-4*
        t1*t34*grada*t35*t59*t58*t49*t54*t39+3.359359811638531E-4*
        t1*t34*grada*t35*t59*t48*t49*t54*t39+3.359359811638531E-4*
        t1*t34*grada*t35*t59*t96*t49*t39-4.479146415518041E-4*t1*t34*
        grada*t35*t97*t58*t49*t39-4.479146415518041E-4*t1*t34*grada*
        t35*t97*t48*t49*t39+3.359359811638531E-4*t1*t34*grada*t35*
        t59*t58*t49*t33*t38*t39+3.359359811638531E-4*t1*t34*grada*
        t35*t59*t48*t49*t33*t38*t39+1.333333333333333*t1*grada*t50*
        t58*t49*t37*t39+1.333333333333333*t1*grada*t50*t48*t49*t37*
        t39+t189+t188+t187+t186+t185+t184+t183+t182+t181+t180+t179+
        t178+t177+t176+t175+t174+t173+t172;
    d3fdrarbgb = -1.0*t1*t37*t39*t49*t6*t96*gradb+1.5*t1*
        t37*t39*t48*t58*t6*t73*gradb-1.0*t1*t37*t39*t49*t54*t58*t6*
        gradb-1.0*t1*t33*t39*t49*t58*t6*gradb-1.0*t1*t37*t39*t48*t49*
        t54*t6*gradb-1.0*t1*t33*t39*t48*t49*t6*gradb-0.5*t1*t111*t37*
        t39*t4*t49*t58*t6-0.5*t1*t111*t37*t39*t4*t48*t49*t6-5.039039717457796E-4*
        t1*t34*gradb*t35*t59*t58*t48*t73*t39+3.359359811638531E-4*
        t1*t34*gradb*t35*t59*t58*t49*t54*t39+3.359359811638531E-4*
        t1*t34*gradb*t35*t59*t48*t49*t54*t39+3.359359811638531E-4*
        t1*t34*gradb*t35*t59*t96*t49*t39-4.479146415518041E-4*t1*t34*
        gradb*t35*t97*t58*t49*t39-4.479146415518041E-4*t1*t34*gradb*
        t35*t97*t48*t49*t39+3.359359811638531E-4*t1*t34*gradb*t35*
        t59*t58*t49*t33*t38*t39+3.359359811638531E-4*t1*t34*gradb*
        t35*t59*t48*t49*t33*t38*t39+1.333333333333333*t1*gradb*t50*
        t58*t49*t37*t39+1.333333333333333*t1*gradb*t50*t48*t49*t37*
        t39+t207+t206+t205+t204+t203+t202+t201+t200+t199+t198+t197+
        t196+t195+t194+t193+t192+t191+t190;
    d3fdrarbab = -1.0*t1*t37*t39*t49*t6*t96+1.5*t1*t37*t39*
        t48*t58*t6*t73-1.0*t1*t37*t39*t49*t54*t58*t6-0.5*t1*t129*t37*
        t39*t4*t49*t58*t6-1.0*t1*t33*t39*t49*t58*t6-1.0*t1*t37*t39*
        t48*t49*t54*t6-0.5*t1*t129*t37*t39*t4*t48*t49*t6-1.0*t1*t33*
        t39*t48*t49*t6-5.039039717457796E-4*t1*t34*t35*t59*t58*t48*
        t73*t39+3.359359811638531E-4*t1*t34*t35*t59*t58*t49*t54*t39+
        3.359359811638531E-4*t1*t34*t35*t59*t48*t49*t54*t39+3.359359811638531E-4*
        t1*t34*t35*t59*t96*t49*t39-4.479146415518041E-4*t1*t34*t35*
        t97*t58*t49*t39-4.479146415518041E-4*t1*t34*t35*t97*t48*t49*
        t39+3.359359811638531E-4*t1*t34*t35*t59*t58*t49*t33*t38*t39+
        3.359359811638531E-4*t1*t34*t35*t59*t48*t49*t33*t38*t39+1.333333333333333*
        t1*t50*t58*t49*t37*t39+1.333333333333333*t1*t50*t48*t49*t37*
        t39+t225+t224+t223+t222+t221+t220+t219+t218+t217+t216+t215+
        t214+t213+t212+t211+t210+t209+t208;
    d3fdragaga = t132+t236+t131+t235+t234+t233+t232+t126+
        t125+t124+t230+t229-2.257059668810414E-7*t1*t22*t2*t137*t48*
        t49*t38*t39+t227+t123+t122+t226+t121+0.00100780794349*t1*t34*
        t2*t102*t59*t48*t49*t39;
    d3fdragagb = t247+t246+t245+t244+t243+t242+t241+t240+
        t239-2.257059668810414E-7*t1*t22*grada*gradb*t137*t48*t49*
        t38*t39+t238+t237+0.00100780794349*t1*t34*grada*gradb*t102*
        t59*t48*t49*t39;
    d3fdrbgagb = t247+t246+t245+t244+t243+t242+t241+t240+
        t239-2.257059668810414E-7*t1*t22*grada*gradb*t137*t58*t49*
        t38*t39+t238+t237+0.00100780794349*t1*t34*grada*gradb*t102*
        t59*t58*t49*t39;
    d3fdragbgb = t132+t256+t131+t255+t254+t253+t252+t126+
        t125+t124+t251+t250-2.257059668810414E-7*t1*t22*t3*t137*t48*
        t49*t38*t39+t249+t123+t122+t248+t121+0.00100780794349*t1*t34*
        t3*t102*t59*t48*t49*t39;
    d3fdrarara = -1.5*t1*t37*t39*t4*t48*t49*t6*t94-1.5*t1*
        t37*t39*t4*t48*t49*t6*t90+2.25*t1*t37*t39*t4*t48*t6*t73*t84-
        1.5*t1*t37*t39*t4*t49*t54*t6*t84-1.5*t1*t33*t39*t4*t49*t6*
        t84-0.5*t1*t37*t39*t4*t49*t6*(1.666666666666667*t46*t7*t80*
        t83+1.666666666666667*t43*t7*t75*t78+0.83333333333333*(t258+
        t157)*t47*t7-0.18518518518519*t156*pow(t46,3.0)*t7+0.83333333333333*
        (t257+t155)*t44*t7-0.18518518518519*t154*pow(t43,3.0)*t7)+
        2.25*t1*t37*t39*t4*t54*t6*t72*t73+2.25*t1*t33*t39*t4*t6*t72*
        t73-1.5*t1*t39*t4*t48*t49*t6*t69-3.0*t1*t33*t39*t4*t48*t49*
        t54*t6-1.875*t1*t151*t37*t39*t4*pow(t48,3.0)*t6+4.0*t1*t4*
        t50*t48*t49*t37*t54*t39-3.0*t1*t4*t50*t72*t73*t37*t39+2.0*
        t1*t4*t50*t84*t49*t37*t39-4.666666666666666*t1*t4*t85*t48*
        t49*t37*t39+4.0*t1*t4*t50*t48*t49*t33*t39+t170+t169+t168+t167+
        t166+t163+t162+t161+t160+t159+t158+t150+t147+t146;
    d3fdrbrbrb = -1.5*t1*t37*t39*t4*t49*t58*t6*t94-1.5*t1*
        t37*t39*t4*t49*t58*t6*t90-0.5*t1*t37*t39*t4*t49*t6*(1.666666666666667*
        t119*t57*t7*t80+1.666666666666667*t117*t56*t7*t75-0.18518518518519*
        t156*pow(t57,3.0)*t7-0.18518518518519*t154*pow(t56,3.0)*t7+
        0.83333333333333*(t257+t157)*t47*t7+0.83333333333333*(t258+
        t155)*t44*t7)+2.25*t1*t120*t37*t39*t4*t58*t6*t73+2.25*t1*t115*
        t37*t39*t4*t54*t6*t73+2.25*t1*t115*t33*t39*t4*t6*t73-1.5*t1*
        t39*t4*t49*t58*t6*t69-1.875*t1*t151*t37*t39*t4*pow(t58,3.0)*
        t6-3.0*t1*t33*t39*t4*t49*t54*t58*t6-1.5*t1*t120*t37*t39*t4*
        t49*t54*t6-1.5*t1*t120*t33*t39*t4*t49*t6+4.0*t1*t4*t50*t58*
        t49*t37*t54*t39-3.0*t1*t4*t50*t115*t73*t37*t39-4.666666666666666*
        t1*t4*t85*t58*t49*t37*t39+2.0*t1*t4*t50*t120*t49*t37*t39+4.0*
        t1*t4*t50*t58*t49*t33*t39+t170+t169+t168+t167+t166+t163+t162+
        t161+t160+t159+t158+t150+t147+t146;
    d3fdrbrbga = 1.5*t1*t115*t37*t39*t6*t73*grada-2.0*t1*
        t37*t39*t49*t54*t58*t6*grada-2.0*t1*t33*t39*t49*t58*t6*grada-
        1.0*t1*t120*t37*t39*t49*t6*grada-1.0*t1*t103*t37*t39*t4*t49*
        t58*t6-5.039039717457796E-4*t1*t34*grada*t35*t59*t115*t73*
        t39+6.718719623277062E-4*t1*t34*grada*t35*t59*t58*t49*t54*
        t39-8.958292831036083E-4*t1*t34*grada*t35*t97*t58*t49*t39+
        3.359359811638531E-4*t1*t34*grada*t35*t59*t120*t49*t39+6.718719623277062E-4*
        t1*t34*grada*t35*t59*t58*t49*t33*t38*t39+2.666666666666667*
        t1*grada*t50*t58*t49*t37*t39+t189+t188+t187+t186+t185+t184+
        t183+t182+t181+t180+t179+t178+t177+t176+t175+t174+t173+t172;
    d3fdrbrbgb = 1.5*
        t1*t115*t37*t39*t6*t73*gradb-2.0*t1*t37*t39*t49*t54*t58*t6*
        gradb-2.0*t1*t33*t39*t49*t58*t6*gradb-1.0*t1*t120*t37*t39*
        t49*t6*gradb-1.0*t1*t111*t37*t39*t4*t49*t58*t6-5.039039717457796E-4*
        t1*t34*gradb*t35*t59*t115*t73*t39+6.718719623277062E-4*t1*
        t34*gradb*t35*t59*t58*t49*t54*t39-8.958292831036083E-4*t1*
        t34*gradb*t35*t97*t58*t49*t39+3.359359811638531E-4*t1*t34*
        gradb*t35*t59*t120*t49*t39+6.718719623277062E-4*t1*t34*gradb*
        t35*t59*t58*t49*t33*t38*t39+2.666666666666667*t1*gradb*t50*
        t58*t49*t37*t39+t207+t206+t205+t204+t203+t202+t201+t200+t199+
        t198+t197+t196+t195+t194+t193+t192+t191+t190;
    d3fdrbgaga = t259;
    d3fdrbgbga = t259;
    d3fdrbgbgb = t132+t256+t131+t255+t254+t253+t252+t126+
        t135+t124+t251+t250-2.257059668810414E-7*t1*t22*t3*t137*t58*
        t49*t38*t39+t249+t123+t122+t248+t134+0.00100780794349*t1*t34*
        t3*t102*t59*t58*t49*t39;
    d3fdgagaga = 1.354235801286248E-6*t1*t22*t260*t263*t137*
        t13*t38*t39+1.354235801286248E-6*t1*t22*grada*t137*t13*t38*
        t39-3.03291021754875E-10*t1*t261*t260*t102*t262*t13*t52*t39-
        0.00604684766095*t1*t34*grada*t102*t59*t13*t39+0.00201561588698*
        t1*t34*t260*t231*t59*t13*t39;
    d3fdgagagb = 1.354235801286248E-6*t1*t22*t2*gradb*t263*
        t137*t13*t38*t39+t141-3.03291021754875E-10*t1*t261*t2*gradb*
        t102*t262*t13*t52*t39+t140+0.00201561588698*t1*t34*t2*gradb*
        t231*t59*t13*t39;
    d3fdgagbgb = 1.354235801286248E-6*t1*t22*grada*t3*t263*
        t137*t13*t38*t39+t139-3.03291021754875E-10*t1*t261*grada*t3*
        t102*t262*t13*t52*t39+t138+0.00201561588698*t1*t34*grada*t3*
        t231*t59*t13*t39;
    d3fdgbgbgb = 1.354235801286248E-6*t1*t22*t264*t263*t137*
        t13*t38*t39+1.354235801286248E-6*t1*t22*gradb*t137*t13*t38*
        t39-3.03291021754875E-10*t1*t261*t264*t102*t262*t13*t52*t39-
        0.00604684766095*t1*t34*gradb*t102*t59*t13*t39+0.00201561588698*
        t1*t34*t264*t231*t59*t13*t39;


    ds->df1000 += factor*dfdra;
    ds->df0100 += factor*dfdrb;
    ds->df0010 += factor*dfdga;
    ds->df0001 += factor*dfdgb;
    ds->df00001+= factor*dfdgab;

    ds->df2000 += factor*d2fdrara;
    ds->df1100 += factor*d2fdrarb;
    ds->df1010 += factor*d2fdraga;
    ds->df1001 += factor*d2fdragb;
    ds->df0200 += factor*d2fdrbrb;
    ds->df10001+= factor*d2fdraab;
    ds->df01001+= factor*d2fdrbab;
    ds->df0020 += factor*d2fdgaga;
    ds->df0011 += factor*d2fdgagb;
    ds->df0002 += factor*d2fdgbgb;
    ds->df0110 += factor*d2fdrbga;
    ds->df0101 += factor*d2fdrbgb;

    ds->df2010 += factor*d3fdraraga;
    ds->df2001 += factor*d3fdraragb;
    ds->df1101 += factor*d3fdrarbgb;
    ds->df11001 += factor*d3fdrarbab;
    ds->df1020 += factor*d3fdragaga;
    ds->df1011 += factor*d3fdragagb;
    ds->df0111 += factor*d3fdrbgagb;
    ds->df1002 += factor*d3fdragbgb;
    ds->df3000 += factor*d3fdrarara;
    ds->df2100 += factor*d3fdrararb;
    ds->df20001 += factor*d3fdraraab;
    ds->df02001 += factor*d3fdrbrbab;
    ds->df1200 += factor*d3fdrarbrb;
    ds->df1110 += factor*d3fdrarbga;
    ds->df0300 += factor*d3fdrbrbrb;
    ds->df0210 += factor*d3fdrbrbga;
    ds->df0201 += factor*d3fdrbrbgb;
    ds->df0120 += factor*d3fdrbgaga;
    ds->df0102 += factor*d3fdrbgbgb;
    ds->df0030 += factor*d3fdgagaga;
    ds->df0021 += factor*d3fdgagagb;
    ds->df0012 += factor*d3fdgagbgb;
    ds->df0003 += factor*d3fdgbgbgb;
}
